本資料は Covid19 Japan が独自に収集している陽性者のデータ(個票データ)を集計・可視化するためのサンプルコード例です。データ収集処理のためのソースは ソースリポジトリ(GitHub)、収集・処理したデータは データリポジトリ(GitHub) にて公開1 ライセンスに関してはリンク先にて確認してください。 されています。
データリポジトリ(GitHub) にて公開されているデータはJSON形式ですので、読み込みには jsonlite パッケージが必要です2 tidyverse パッケージと共にインストールされます。インストールされていない場合は以下のコードでインストールしてください。。 install.packages(“jsonlite”)
陽性と判定された陽性判定者(以降、陽性者と記述)単位で記録されている個票データは データリポジトリ(GitHub) から直接読み込めます3 時間帯によりデータが揃っていない場合があります。その場合は、時間帯を変更するか、前日までのデータを利用してください。 が、必ず raw.githubusercontent.com のパスを使ってください。github.com パスでは読み込めませんので注意してください。 コードは read.csv や readr::read_csv を使う場合と同様で下記のようになります。
df <- "https://raw.githubusercontent.com/reustle/covid19japan-data/master/docs/patient_data/latest.json" %>%
jsonlite::fromJSON()
読み込んだ個票データを skimr パッケージを利用して要約すると下記のようになります4 数値と論理型を除き文字型として扱われますので適切な型に変換する必要があります。。
df %>%
skimr::skim()
Data summary
| Name | Piped data |
| Number of rows | 445747 |
| Number of columns | 17 |
| _______________________ | |
| Column type frequency: | |
| character | 15 |
| logical | 1 |
| numeric | 1 |
| ________________________ | |
| Group variables | None |
Variable type: character
| skim_variable | n_missing | complete_rate | min | max | empty | n_unique | whitespace |
|---|---|---|---|---|---|---|---|
| patientId | 0 | 1.00 | 1 | 16 | 0 | 437635 | 0 |
| dateAnnounced | 0 | 1.00 | 10 | 10 | 0 | 403 | 0 |
| gender | 335027 | 0.25 | 1 | 1 | 0 | 2 | 0 |
| detectedPrefecture | 0 | 1.00 | 3 | 15 | 0 | 49 | 0 |
| patientStatus | 435380 | 0.02 | 8 | 23 | 0 | 8 | 0 |
| mhlwPatientNumber | 445298 | 0.00 | 1 | 11 | 0 | 434 | 0 |
| prefecturePatientNumber | 328302 | 0.26 | 5 | 20 | 0 | 117436 | 0 |
| residence | 343117 | 0.23 | 1 | 38 | 0 | 1429 | 0 |
| relatedPatients | 433160 | 0.03 | 2 | 259 | 0 | 7470 | 0 |
| knownCluster | 443213 | 0.01 | 3 | 88 | 0 | 235 | 0 |
| detectedCityTown | 416714 | 0.07 | 2 | 22 | 0 | 667 | 0 |
| cityPrefectureNumber | 417015 | 0.06 | 1 | 34 | 0 | 28723 | 2 |
| deceasedDate | 437541 | 0.02 | 10 | 10 | 0 | 353 | 0 |
| deceasedReportedDate | 444462 | 0.00 | 10 | 62 | 0 | 209 | 0 |
| deathSourceURL | 444605 | 0.00 | 14 | 123 | 0 | 659 | 0 |
Variable type: logical
| skim_variable | n_missing | complete_rate | mean | count |
|---|---|---|---|---|
| confirmedPatient | 0 | 1 | 0.98 | TRU: 437634, FAL: 8113 |
Variable type: numeric
| skim_variable | n_missing | complete_rate | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|
| ageBracket | 335140 | 0.25 | 38.01 | 20.69 | 0 | 20 | 30 | 50 | 100 | ▇▇▅▂▁ |
Subset of JSON data
| patientId | dateAnnounced | ageBracket | gender | detectedPrefecture | |
|---|---|---|---|---|---|
| 1 | 15 | 2020-01-15 | 30 | M | Kanagawa |
| 2 | TOK1 | 2020-01-24 | 40 | M | Tokyo |
| 3 | TOK2 | 2020-01-25 | 30 | F | Tokyo |
| 4 | 18 | 2020-01-26 | 40 | M | Aichi |
| … | NA | NA | … | NA | NA |
| 445744 | -1 | 2021-03-06 | NA | NA | Yamaguchi |
| 445745 | -1 | 2021-03-06 | NA | NA | Fukushima |
| 445746 | -1 | 2021-03-06 | NA | NA | Fukushima |
| 445747 | -1 | 2021-03-06 | NA | NA | Fukushima |
個票の データフォーマット(GitHub) を参考に適切な型5 日付は日付型、属性は因子型に変換しておくと処理しやすくなります。 に変換するとともに不要な項目を外しておきます。
また、記録には重複や未確定なデータが含まれていますので、これらを外しておく6 confirmedPatient == TRUE でフィルタリングします。 ことも必要です。
dft <- df %>%
dplyr::select(
patientId, date = dateAnnounced, gender, pref = detectedPrefecture,
patientStatus, knownCluster, confirmedPatient, ageBracket
) %>%
dplyr::filter(confirmedPatient == TRUE) %>%
dplyr::mutate(
date = lubridate::as_date(date),
gender = forcats::as_factor(gender),
pref = stringr::str_to_lower(pref),
patientStatus = forcats::as_factor(patientStatus),
cluster = dplyr::if_else(!is.na(knownCluster), TRUE, FALSE),
ageBracket = forcats::as_factor(ageBracket)
)
変換結果を要約すると下記のように適切に変換されれていることが分かります。
| Name | Piped data |
| Number of rows | 437634 |
| Number of columns | 9 |
| _______________________ | |
| Column type frequency: | |
| character | 3 |
| Date | 1 |
| factor | 3 |
| logical | 2 |
| ________________________ | |
| Group variables | None |
Variable type: character
| skim_variable | n_missing | complete_rate | min | max | empty | n_unique | whitespace |
|---|---|---|---|---|---|---|---|
| patientId | 0 | 1.00 | 1 | 16 | 0 | 437634 | 0 |
| pref | 0 | 1.00 | 3 | 15 | 0 | 49 | 0 |
| knownCluster | 435129 | 0.01 | 3 | 88 | 0 | 233 | 0 |
Variable type: Date
| skim_variable | n_missing | complete_rate | min | max | median | n_unique |
|---|---|---|---|---|---|---|
| date | 0 | 1 | 2020-01-15 | 2021-03-06 | 2020-12-27 | 403 |
Variable type: factor
| skim_variable | n_missing | complete_rate | ordered | n_unique | top_counts |
|---|---|---|---|---|---|
| gender | 327962 | 0.25 | FALSE | 2 | M: 61370, F: 48302 |
| patientStatus | 435100 | 0.01 | FALSE | 8 | Hos: 1261, Dec: 372, Hom: 315, Dis: 283 |
| ageBracket | 328059 | 0.25 | FALSE | 12 | 20: 29433, 30: 19042, 40: 16089, 50: 14139 |
Variable type: logical
| skim_variable | n_missing | complete_rate | mean | count |
|---|---|---|---|---|
| confirmedPatient | 0 | 1 | 1.00 | TRU: 437634 |
| cluster | 0 | 1 | 0.01 | FAL: 435129, TRU: 2505 |
個票データの都道府県(detectedPrefecture)を利用して都道府県に関するデータ7 著者が作成し Gist で公開しているもの。 をマージすることで、地方などの地域を切り口とした集計が可能になります。 個票データを除く他のデータ8 Google感染予測 など とのマージも考慮して項目名や型を適宜変換しておきます。
prefs <- "https://gist.githubusercontent.com/k-metrics/9f3fc18e042850ff24ad9676ac34764b/raw/f4ea87f429e1ca28627feff94b67c8b2432aee59/pref_utf8.csv" %>%
readr::read_csv() %>%
dplyr::mutate(
# Googleの予測データと結合を考慮してコード体系を合わせておく
japan_prefecture_code = paste0("JP-", `コード`)
) %>%
dplyr::select(
# Googleの予測データと結合するために名称を変更する
japan_prefecture_code, prefecture_name = pref,
# 日本語の変数名は扱いにくいので英語名に変更する
pref = `都道府県`, region = `八地方区分`, pops = `推計人口`
) %>%
dplyr::mutate(
# 水準ごとに表示させるために因子化する(あらかじめデータをコード順に
# 並べておくことが因子化の際のポイントのひとつ)
japan_prefecture_code = forcats::fct_inorder(japan_prefecture_code),
pref = forcats::fct_inorder(pref),
region = forcats::fct_inorder(region),
pops = as.integer(pops)
)
pops)は千人単位になっています。
前出の個票データと都道府県データを結合します。結合の際に因子型変数の水準が欠損値となっているものは適宜 forcats パッケージを用いて水準を割り当てておきます。
x <- dft %>%
dplyr::left_join(prefs, by = c("pref" = "prefecture_name")) %>%
dplyr::select(-pref) %>%
dplyr::rename(pref = pref.y) %>%
# 因子型の欠損値を水準化しておく
dplyr::mutate(
japan_prefecture_code = forcats::fct_explicit_na(japan_prefecture_code,
na_level = "JP-48"),
pref = forcats::fct_explicit_na(pref, na_level = "空港検疫"),
region = forcats::fct_explicit_na(region, na_level = "空港検疫"),
gender = forcats::fct_explicit_na(gender, na_level = "非公表"),
ageBracket = forcats::fct_explicit_na(ageBracket, na_level = "非公表"),
patientStatus = forcats::fct_explicit_na(patientStatus,
na_level = "Unknown")
)
完成した個票データを用いて様々な集計を行います。
日付(date)を用いた単純な日次集計を行います。この集計結果は 厚生労働省オープンデータ の陽性者数データに相当9 厚生労働省の場合、途中から個票からの集計を取りやめ各自治体が公表している集計値を積上げたものを公開しています。 します。
集計には dplyr::group_by と dplyr::summarize(n = n())10 dplyr::count を用いても同様の結果を得ることができます。 を用いることで簡単に計数することができます。ただし、個票データが存在しない日付は計数ができませんので暗黙の欠落(欠損)となり明示的な欠落(欠損)に変換する11 Turns implicit missing values into explicit missing values. 必要があります。この処理は tidyr::complete を用いることで簡単にできます。
japan_daily <- x %>%
dplyr::group_by(date) %>%
dplyr::summarise(n = dplyr::n()) %>%
dplyr::ungroup() %>%
tidyr::complete(
date = seq.Date(from = min(date), to = max(date), by = "day"),
fill = list(n = 0L)
)
dplyr::group_by と dplyr::summarize を行った後に続けて処理を行う場合には dplyr::ungroup で必ずアングループしてください。
前日差(diff)・累計(cumsum)・移動平均(7日(ma7)ならびに28日(ma28))12 lagdiff, ma7, ma28 は以下のような関数として定義してあるます。 も求めます。
lagdiff <- function(n) { n - dplyr::lag(n, default = 0L) }
ma7 <- function(n) { zoo::rollmeanr(n, k = 7L, na.pad = TRUE) }
ma28 <- function(n) { zoo::rollmeanr(n, k = 28L, na.pad = TRUE) }
japan_daily <- japan_daily %>%
dplyr::mutate(
diff = lagdiff(n), # 前日差
cum = cumsum(n), # 累計
ma7 = ma7(n), # 移動平均(7日)
ma28 = ma28(n) # 移動平均(28日)
)
都道府県ごとの日次集計、年代別の日次集計など複数の変数の水準ごとの集計をクロス集計と呼びます。このクロス集計も単純集計と同様に dplyr::group_by と dplyr::summrize(n = n())13 単純集計と同様に dplyr::count を用いても同様の結果を得ることができます。 で計数することができます。
地方区分(region)を dplyr::group_by に追加指定するだけでクロス集計ができます。暗黙の欠落(欠損)を明示的な欠落(欠損)に変換する場合も tidyr::complte に地方区分(region)を追加するだけ済みます。
region_daily <- x %>%
dplyr::group_by(date, region) %>%
dplyr::summarise(n = dplyr::n()) %>%
dplyr::ungroup() %>%
tidyr::complete(
date = seq.Date(from = min(date), to = max(date), by = "day"), region,
fill = list(n = 0L)
) %>%
dplyr::group_by(region) %>%
dplyr::mutate(
diff = lagdiff(n),
cum = cumsum(n),
ma7 = ma7(n),
ma28 = ma28(n)
) %>%
dplyr::ungroup()
集計結果を可視化してみます。
日次集計データ(japan_daily)の単日陽性者数(n)・累計陽性者数(cum)・単日移動平均(7日)(ma7)を描きます。単日陽性者数は棒グラフ、累計陽性者数と単日移動平均は折線グラフで描きますが、累計の最大値は単日や移動平均と2桁異なりますので二軸のグラフとして描く必要があります。
title <- "Tufte's Style"
xlab <- ""
ylab <- ""
sec_scale <- 50 # 縦二軸用のスケーリング値
dbreaks <- "2 month" # 横軸目盛りのラベル間隔
dlabels <- "%y-%m" # 横軸目盛りラベルのフォーマット
axis1 <- "陽性者数(灰)・同移動平均(濃灰)"
axis2 <- "累計陽性者数(濃緑)"
japan_daily %>%
ggplot2::ggplot(ggplot2::aes(x = date)) +
ggplot2::geom_bar(ggplot2::aes(y = n), stat = "identity", width = 1.0,
fill = "dark gray", alpha = 0.75) +
ggplot2::geom_line(ggplot2::aes(y = ma7), linetype = "solid",
colour = "gray10", size = 0.35) +
# 第二軸を利用するグラフを描画する際はスケーリング調整する
ggplot2::geom_line(ggplot2::aes(y = cum / sec_scale),
colour = "dark green", size = 0.75) +
# 横軸表示の指定
ggplot2::scale_x_date(date_breaks = "1 month", date_labels = "%y/%m") +
# 二軸表示のための軸属性の指定
ggplot2::scale_y_continuous(
# 第一軸のラベル(スケールは自動調整)
name = axis1,
# 第二軸の指定(第一軸にスケーリング値をかけたスケール)
sec.axis = ggplot2::sec_axis(~ . * sec_scale,
name = axis2)) +
ggthemes::theme_tufte() +
ggplot2::scale_x_date(date_breaks = dbreaks, date_labels = dlabels) +
ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
x = xlab, y = ylab)
【全国】陽性者数の推移(単日/累計/移動平均(7日))
Draw with default theme.
上グラフは
ggthemes::theme_tufteを適用し Tufte スタイルで描画しています。右側のデフォルトテーマ(ggplot2::theme_gray)で描画したグラフと比べるとデータ・インクレシオが高いことが分かります。
ggthemes::theme_tufte は凡例を非表示にできませんので ggplot2::theme を利用してTufteスタイルに近い設定にします。
title <- "地方別陽性者数の推移"
xlab <- ""
ylab <- ""
sec_scale <- 50 # 縦二軸用のスケーリング値
ncol <- 3
dbreaks <- "3 month"
dlabels <- "%y-%m"
axis1 <- "陽性者数(棒)・同移動平均(細線)"
axis2 <- "累計陽性者数(太線)"
region_daily %>%
dplyr::mutate(key = region) %>%
ggplot2::ggplot(ggplot2::aes(x = date)) +
ggplot2::geom_bar(ggplot2::aes(y = n, fill = key), stat = "identity",
alpha = 0.25, width = 1.0) +
ggplot2::geom_line(ggplot2::aes(y = ma7, colour = key),
linetype = "solid", size = 0.25) +
ggplot2::geom_line(ggplot2::aes(y = cum / sec_scale, colour = key)) +
ggplot2::scale_x_date(date_breaks = dbreaks, date_labels = dlabels) +
ggplot2::theme_gray() +
ggplot2::theme(panel.grid = ggplot2::element_blank(),
panel.background = ggplot2::element_blank(),
strip.background = ggplot2::element_blank(),
legend.position = "none") +
ggplot2::facet_wrap(~ key, ncol = ncol, scales = "free_y") +
ggplot2::scale_y_continuous(
name = axis1,
sec.axis = ggplot2::sec_axis(~ . * sec_scale, name = axis2)) +
ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
x = xlab, y = ylab)
title <- "都道府県別陽性者数の推移"
xlab <- ""
ylab <- ""
sec_scale <- 50 # 縦二軸用のスケーリング値
ncol <- 4
dbreaks <- "3 month"
dlabels <- "%y-%m"
axis1 <- "陽性者数(棒)・同移動平均(細線)"
axis2 <- "累計陽性者数(太線)"
pref_daily %>%
dplyr::mutate(key = pref) %>%
ggplot2::ggplot(ggplot2::aes(x = date)) +
ggplot2::geom_bar(ggplot2::aes(y = n, fill = key), stat = "identity",
alpha = 0.25, width = 1.0) +
ggplot2::geom_line(ggplot2::aes(y = ma7, colour = key),
linetype = "solid", size = 0.25) +
ggplot2::geom_line(ggplot2::aes(y = cum / sec_scale, colour = key)) +
ggplot2::scale_x_date(date_breaks = dbreaks, date_labels = dlabels) +
ggplot2::theme_gray() +
ggplot2::theme(panel.grid = ggplot2::element_blank(),
panel.background = ggplot2::element_blank(),
strip.background = ggplot2::element_blank(),
legend.position = "none") +
ggplot2::facet_wrap(~ key, ncol = ncol, scales = "free_y") +
ggplot2::scale_y_continuous(
name = axis1,
sec.axis = ggplot2::sec_axis(~ . * sec_scale, name = axis2)) +
ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
x = xlab, y = ylab)